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INTRODUCTION 


This report is the second annual summary of the research program of the 
Geodynamics Branch. The branch is located within the Laboratory for Earth 
Sciences of the Applications Directorate of the Goddard Space Flight Center. 
The research activities of the branch staff cover a broad spectrum of geoscience 
disciplines including space geodesy, geopotential field modeling, tectono- 
physics, and dynamic oceanography. The NASA programs which are supported by 
the work described in this document, include the Geodynamics and Ocean Programs, 
the Crustal Dynamics Project, the proposed Ocean Topography Experiment (TOPEX), 
and the Geopotential Research Mission (GRM). The reports highlight the inves- 
tigations conducted by the Geodynamics Branch Staff during calendar year 1983. 
The individual papers are grouped into chapters on Crustal Movements, Global 
Earth' Dynamics, Gravity Field Model Development, Sea Surface Topography, and 
Advanced Studies. Further information on the activities of the branch or the 
particular research efforts described here can be obtained through the branch 
office or from the individual staff members. 
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CHAPTER 1 
CRUSTAL MOVEMENTS 


OVERVIEW 


Studies of crustal movements are devoted to the interrelated goals 
of determining tectonic plate motions, understanding their driving 
forces and the nature of plate interactions, and predicting the 
occurrence of earthquakes. The analysis of global and regional 
scale crustal movements provide information on seafloor spreading 
and continental drift, the creation and destruction of crustal 
material, and the structure, rheology, and dynamics of the earth's 
interior. Crustal deformation in seismically active zones provides 
information on strain accumulation and release, the state of stress 
in the earth's crust, and the mechanisms responsible for earth- 
quakes. The work of the Geodynamics Branch includes activities in 
both the Crustal Dynamic Project which involves space geodetic 
observations and in the Geodynamics Research Program which involves 
basic studies into the dynamic processes of the earth. 

This section of the Geodynamics Branch Report includes papers on 
the Measurement of Crustal Movements, Modeling of Earthquake Related 
Crustal Deformations, Modeling of Sublithospheric Stress Fields, 
and the Development of Satellite Tracking Facilities for Space 
Geodetic Measurements. The six reports are summarized as follows: 
1) The first report is devoted to discussing recent results obtained 
from the San Andreas Fault Experiment. This decade long experiment 
has utilized satellite laser ranging techniques to monitor the 
broad motion of the tectonic plates comprising the San Andreas 
Fault System. 2) The second report focuses on the modeling of 
crustal deformation subsequent to an earthquake. The analyses are 
based on the physics of stress redistribution with time and distance 
from the seismic event and yields information on the rheology and 
structure of the earth's interior. 3) The third report discusses 
the results of Satellite Laser Ranging System collocation experi- 
ments performed to assess the laser system's performance prior to 
field deployment of the laser or during periods when the laser 
system's performance is degraded. 4) The fourth report describes 
how two different measurement techniques, namely Very Long Baseline 
Interferometry, and Satellite Laser Ranging compare in the measure- 
ment of distances between several locations in the continental U.S. 
5) In the fifth report, the sublithospheric stress field in central 
Europe is derived from an analysis of the satellite-determined 
gravity field. The computed stress field is interpreted in terms 
of the regional seismicity, and geologic structure, and the pheno- 
menon of gravity sliding of seismotectonic blocks. 6) The sixth 
report describes the result of error analyses undertaken to assess 
the quality of station coordinates, polar motion and earth rotation 
obtained from the analysis of laser ranging data from Lageos from 
March (1976) to the end of 1982. 

Contributors to this chapter are: D.C. Christodoul idis, S.C. Cohen, 
R. Kolenkiewicz, H.S. Liu, J.W. Ryan, and D.E. Smith. 
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UPDATE: 

SAN ANDREAS FAULT EXPERIMENT 


Demosthenes C. Christodoul idis 
David E. Smith 


OBJECTIVE 

A decade long experiment has utilized satellite laser ranging 
techniques to monitor the broad motion of the tectonic plates 
comprising the San Andreas Fault System. 


BACKGROUND 

Commencing in 1972, mobile laser stations have occupied sites on 
opposite sides and distances from the San Andreas Fault. This San 
Andreas Fault Experiment, (SAFE), has progressed through the up- 
grades made to laser system hardware and an improvement in the 
modeling capabilities of the spaceborne laser targets. Of special 
note is the launch of the LAGEOS spacecraft, NASA's only completely 
dedicated laser satellite in 1976. The results of plate motion 
projected into this 896 km measured line over the past eleven 
years are summarized and intercompared below. 


RECENT ACCOMPLISHMENT 

In satellite laser ranging the observed quantities are the nearly 
instantaneous ranges to a spacecraft and the times at which these 
distances are measured. If the evolution of the orbiting target 
and the rotational variations of the earth are accurately known, 
then the position of the observing station can be calculated in the 
geocentric reference frame used to describe the satellite motion. 

With the launch of LAGEOS at high altitude, a laser target became 
available which was far less sensitive to the poorly known pertur- 
bing forces of short wavelength gravity and atmospheric drag. 
This is in marked contrast to the modeling problems found with 
BE-C, the spacecraft used for SAFE prior to launch of LAGEOS. 
With the development of a worldwide network of more advanced 
laser systems, the trajectory accuracies achievable for LAGEOS 
approach the quality of the laser tracking itself. Improved long 
wavelength gravity fields and highly accurate parameters to orient 
the earth's pole and rotation evolved from the analysis of the 
LAGEOS data. Thus, the requirement for specially "tailored" SAFE 
solutions as was done using BE-C data was no longer a necessity. 
The results obtained both from BE-C and Lageos data analyses are 
quite consistent and of a great deal of geophysical interest. 
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SIGNIFICANCE 

A summary of the changes in the SAFE baseline lengths are shown in 
Figure 1. This figure is made somewhat complex due to the fact 
that the southern site for the experiment was changed from San 
Diego to Monument Peak, a distance of 50 km, during the course -of 
the Crustal Dynamic's Project. An adequate survey tie between San 
Diego and Monument Peak, comparable to the accuracy of the SLR 
measurements is lacking. Therefore, we have not attempted to 
present the length changes projected onto a single line, but rather, 
show changes in length over consecutive measures of the actual 
lines which were observed. Two independent sets of results from 
LAGEOS are compared to those from BE-C. One set of measurements 
(from 1976 to 1979) is coincident between the two satellites and 
directly compared. 

Baseline lengths (referenced to an arbitrary set of origins) which 
were obtained from the laser tracking data are shown in Figure 1. 
A weighted fit of a straight line to each set of values was calcu- 
lated to assess the rate of change in the length of the baselines. 
Since all these lines address the same broad tectonic motions, 
they can be combined to monitor an eleven year history of tectonic 
motion along the San Andreas Fault. 

If the independent SLR rate measurements ere assumed to be measures 
of the same constant tectonic rate, the weighted average rate for 
the last eleven years is -6.4 ± 1.0 cm/year with all measures 
agreeing to within their uncertainty. There is also a high level 
of agreement for the rates obtained using BE-C with those gotten 
from Lageos. 

The satellite laser results can be compared to an independent 
estimate of tectonic motion. Minster and Jordan (1978) have made a 
global model of plate motion which reflects average plate motion 
(assuming rigid plates) over millions of years. Their model incor- 
porates numerous data sources including earthquake slip vectors, 
transform fault directions and the spreading rate deduced from 
magnetic reversals at oceanic ridges. When using this model, the 
estimated rate of change for the San Diego to Quincy baseline is 
expected to be -5.6 cm/year. This value is within the formal error 
of our estimate of -6.4 ± 1.0 cm/year although it is hard to ascribe 
an uncertainty directly to the Minster and Jordan estimate. 


FUTURE EMPHASIS 

Satellite laser ranging has reached a new maturity with the launch 
of LAGEOS and the deployment of a globally distributed network of 
third generation trackers. The analysis of laser data has produced 
an eleven year history of inter-plate motion along the San Andreas 
Fault in California. The broad scale rate of motion between the 
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North American and Pacific Plates projected in the direction of 
Quincy to San Diego shows a baseline shortening of 6.4 ± 1.0 
cm/years. This rate is in good agreement with the value predicted 
by the global tectonic model of Minster and Jordan (1978). Analysis 
techniques to improve the temporal resolution of the SLR measures 
are under development with preliminary results indicating that our 
standard annual solutions can be reduced to ' monthly or shorter 
samples of this inter-plate, motion. 
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Figure 1, San Andreas Fault Experiment 
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CRUSTAL DEFORMATION AND EARTHQUAKES 


Steven C. Cohen 


OBJECTIVE 

The purpose of this work is to develop an understanding of how the 
earth's surface deforms during the cycleof stress accumulation and 
release along major faults. 


BACKGROUND 

The spatial and temporal patterns of crustal deformation that 
accompany various phases of the earthquake cycle depend on the 
dynamical and structural properties of the stressed material in 
fault zones. The presence of viscoelastic material below the 
seismogenic zone in either the asthenosphere or in an intractrustal 
low viscosity zone has a pronounced effect on surface movement. 
The occurrence of an earthquake transfers stress from near surface 
elastic brittle material to deeper semi-fluid material. Viscoe- 
lastic flow of subsurface material in" response to the stressing 
produces surface deformation due to the mechanical coupling between 
layers. The specific pattern that develops depends on the fault 
zone geometry, the thickness of the lithosphere and asthenosphere, 
and the viscosity of the semi-fluid material. Several different 
viscoelastic models of the earthquake cycle have been proposed for 
strike-slip geometries. These models, shown in Figure 1, differ in 
whether the asthenosphere is taken to be thin or thick, and whether 
there are lateral variations in the mechanical properties of 
subsurface rock. Not surprisingly, the predicted deformation 
patterns associated with these models have significant dissimilari -. 
ties. However the model differences arise not only from differ- 
ences in physics but also from mathematical treatments employed to 
reduce the governing differential equations to a tractable form. 
The introduction of mathematical approximations can be largely 
overcome by employing finite element techniques uniformly to all 
the models. Numerical simulation of the earthquake cycle using 
the finite element approach permits model resolution and thereby 
improves our understanding of the details of the crustal deforma- 
tion process. 


RECENT ACCOMPLISHMENTS 

It has been demonstrated that the surface crustal deformation pat- 
tern is sensitive mostly to the rheology of the material that lies 
below the seismic slip plane in a volume whose radius is a few 
times the fault depth. If this material is viscoelastic, 
the surface deformation pattern resembles that of an elastic 
layer over a viscoelastic half-space. However when the thickness 
or breadth of the viscoelastic material is comparable in size to 
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the fault depth, then the surface deformation pattern is altered 
by the geometry. In this case, geodetic measurements are poten- 
tially useful for studying the subsurface structure. Model dif- 
ferences are most detectable at intermediate distances from the 
fault {10^ - 10^ km) where the bFoad-sc-ale features of viscoelas- 
tic flow have their clearest resolution. These differences are 
most pronounced shortly after an earthquake when the viscoelastic 
mechanisms generate signals much larger than those due to elastic 
mechanisms acting alone. 

In an investigation of the crustal deformation associated with a 
thin channel asthenosphere it has been found that displacements are 
reduced from those computed for a half-space asthenosphere (see 
Figure 2). This result demonstrates that a previous finding by 
other workers that displacements are enhanced when flow is confined 
to a thin channel is based on several invalid approximations. The 
major predictions of the finite-element model are that the near- 
field postseismic displacements and strain rates are less than 
those for a half-space asthenosphere and that the postseismic strain 
rates at intermediate distances are greater (in -magnitude).” However 
the finite width of the asthenosphere ceases to have a significant 
impact on the crustal deformation pattern when i ts. magnitude exceeds 
about three'lithosphere thicknesses. 


FUTURE EMPHASiS 

The theoretical modeling aspects of the current research are now 
complete. The main thrust of current and future work is to apply 
the models to observed crutal movements and to extended finite 
element techniques to other crustal -deformation processes such as 
continental collisions, postglacial rebound, and plate bending. 
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Figure 1. Models Used in Earthquake Cycle Simulations: 

a. Elastic Half-space Mode! 

b. Viscoelastic Half-space Asthenosphere Model 

c. Viscoelastic Thin-channel Asthenosphere Model 

d. Viscoelastic Varying Lithosphere Mode! 

e. Viscoelastic Inclusion Model 
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Figure 2. Comparison of Patterns of Displacement Versus Distance frpm the Fault 
for Thin Channel Analytic Model and Finite Element Model. 
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RESULTS OF LASER RANGING COLLOCATIONS 
DURING 1983' 

Ronald Kolenkiewicz 


OBJECTIVE 

The objective of laser ranging collocations is to compare the 
ability of two satellite laser ranging systems, located in the 
vicinity of one another, to measure the distance to an artificial 
earth satellite in orbit over the sites. The similar measurement 
of this distance is essential before a new or modified laser system 
is deployed to worldwide locations in order to gather the data neces- 
sary to meet the scientific goals of the Crustal Dynamics Project. 


BACKGROUND 

For the direct study of crustal movements, position accuracies of a 
few centimeters are required over distances of hundreds or thousands 
of kilometers. Classical methods are inherently inadequate for 
this purpose. The space approach is to use an extraterrestrial 
reference system for determining the position of points on the 
Earth's surface. 

One of the techniques applied by the Crustal Dynamics Project to do 
this uses laser ranging to satellites (i.e., measuring the time of 
flight of very short laser pulses to a retroreflector above the 
Earth's surface). This requires a high-altitude satellite whose 
orbital position can be determined to appropriate accuracy. The 
Laser Geodynamic Satellite (Lageos) is an example of a retrore- 
flector equipped satellite in orbit, which is high enough to be 
free from serious perturbations caused by effects that are diffi- 
cult to model (i.e., the complex gravity field structure attribut- 
able to density inhomogeneities in -the Earth and variations in 
atmospheric drag). By knowing the satellite ephemeris, the posi- 
tion of a laser ranging station can be determined in an earth 
center-of-mass coordinate system. 

In order to be certain the laser systems are operating properly, 
they are periodically compared with each other. These comparisons 
or collocations are performed by locating the lasers side by side 
when they track the same satellite during the same time or pass. 
The data is then compared to make sure the lasers are giving essen- 
tially the same range results. Results of the three collocations 
performed during 1983 are given. 
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RECENT ACCOMPLISHMENTS AND SIGNIFICANCE 


During 1983 three collocation experiments were conducted viz: 

1. Mobile laser system 4 (M08LAS 4) was compared with mobile 

laser system 7 (MOBLAS 7) at Greenbelt, MD in May and June. 

The purpose of this collocation was to check out the newly 

upgraded MOBLAS 4 before deployment to Monument Peak, CA. 

2. Mobile laser system 8 (MOBLAS 8) was compared with transport- 
able laser system (TLRS 1) at Quincy, CA in July and August. 
This was a routine test to see how both systems were per- 
formi ng. 

3. MOBLAS 4 was compared with TLRS 1 at Monument Peak, CA in 

October and November. MOBLAS 4 had been struck by lightning 
and its performance was being evaluated before it resumed its 
tracking operations. 

A figure showing the mean residual difference (in centimeters), per 
pass of the satellite, between the two tracking systems is shown 
for each of the collocations (Figures 1,2, and 3). 

In all cases the satellite tracked was Lageos. For example the 

collocation between MOBLAS 7 and MOBLAS 4 at Greenbelt, MD consisted 
of tracking the Lageos satellite by both lasers on thirteen separate 
passes. During the first pass. May 25, 1983, the MOBLAS 7 laser 
obtained 1050 range data points and MOBLAS 4 got 1395 points. Using 
these data an orbit was determined. The mean of the residuals 
(observed minus computed) about this orbit for MOBLAS 7 was -1.46 
cm and for MOBLAS 4 was 0,84 cm. The difference of these means 
(MOBLAS 7 minus MOBLAS 4)' is -2.30. In other words, for this 
pass, MOBLAS 7 was measuring shorter than MOBLAS 4 by 2.30 centi- 
meters. The mean of these differences for the entire thirteen 
pass collocation was -2.95 4.30 cm. 

At this stage of evolution of the laser systems a criterion for a 
successful collocation requires that the mean of the differences 
obtained during the entire collocation be less then 3 centimeters. 
MOBLAS 4 marginally passed collocation and was transported to 
Monument Peak to join the satellite laser tracking network. 

Figure 2 presents the Quincy collocation which yielded an RMS mean 
residual difference for eighteen passes of Lageos of -1.88 +_ 3.29 
cm. MOBLAS 8 and TLRS 1 are essentially yielding the same values 
ranging to the satellite and are considered to have a successful 
col location. 
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Figure 3 shows the Monument Peak collocation for ten passes of 
Lageos. The RMS mean difference for these passes is 4.06 _+ 1.60 
cm. This is not a successful collocation between these two systems 
since there is a rather large bias (MOBLAS 4 longer by 4 cm) in 
the tracking data. This bias could be due to a number of things. 
Among them: (1) Either laser could.be malfunctioning; (2) a survey 

error to the calibration board or between the lasers could have 
been made. MOBLAS 4 continues to track at Monument Peak and TLRS 1 
has been deployed to other locations. Although useful systems they 
are sti.ll suspect until further analysis of the data has been 
completed. 


FUTURE EMPHASIS 

An attempt will be made to get better agreement between collocated 
laser ranging systems. This will be done by engineering changes to 
upgrade the systems, tighter constraints on the local surveys, and a 
review of current calibration and operation procedures. 


ACKNOWLEDGMENTS 
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Figure 1 . Greenbelt Collocation MOBLAS 7 Minus MOBLAS 4 
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Figure 2. Quincy Collocation MOBLAS 8 Minus TLRS 1 
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Figure 3. Monument Peak Collocation MOBLAS 4 Minus TLRS 1 
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A COMPARISON BETWEEN LAGEOS LASER RANGING 

AND 

VLB I DETERMINED BASELINES 


Ronald Kolenkiewicz 
Janies W. Ryan 


OBJECTIVE 

The objective of this research is to determine how well two 
independent measurement techniques, Lageos satellite laser ranging 
(SLR), and very long baseline interferometry (VLBI) compare in the 
measurement of distances (or baselines) between several locations 
in the continental U.S. 


BACKGROUND 

The NASA Crustal Dynamics Project applies space methods and techno- 
logy to advance the scientific understanding of earth dynamics. 
An important item in this understanding involves precise determina- 
tion of position or rate of change of position in three dimensions 
and at places separated by hundreds-to-thousands of kilometers on 
the surface of the earth. The rates of change of position expec- 
ted from geophysical and seismological considerations are approxi- 
mately several centimeters per year, implying that over reasonable 
time scales, observational systems must be able to detect position 
changes of a few centimeters. These are the measurements that the 
Crustal Dynamics Project plans to collect, analyze, and interpret. 
Currently only two systems are capable of making these measurements, 
SLR and VLBI. In order to insure that the measurements are correct 
it is desirable to compare the results obtained by these systems. 

Since October 1979 the SLR and VLBI elements of NASA's Crustal 
Dynamics Project have conducted coordinated experiments in which 
each technique independently measures the same baseline lengths. 
Sites occupied by the systems include: Westford, MA; Fort Davis, 
TX; Platteville, CO; and Owens Valley, Goldstone, Monument Peak, 
Quincy, and Pasadena, CA. Of the possible 28 baselines between 
these sites, 22 have been measured by both systems. 


RECENT ACCOMPLISHMENTS AND RESULTS 

Analysis of the data to produce the 22 sets of baselines for 
comparison has been completed. For SLR, Lageos laser data was 
analyzed in 92 30-day segments (ARCS) utilizing tracking data from 
a globally-distributed network and did not require simultaneous 
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tracking between observing stations. The parameters solved for in 
these arcs are the Lageos orbit, an along-track satellite acceleration, 
a solar radiation pressure coefficient, and polar motion and A1-UT1 
every 5 days fixing the first A1-UT1 value in each arc. A value of 
GM and the latitude, longitude and height of each of the tracking 
sites (except the latitude and longitude of .Arequipa-i Peru and the 
latitude of 'Hawaii ) were adjusted in common from all data. The 
force models employed included GEM-L2 gravity model with a recovered 
GM value of 398600.44 km^/sec^, Wahr's model for solid-earth tides 
and the Schwiderski Model was used to model the ocean tides. 

The VLBI measurements were made during 24 days since August 1981. The 
observing sessions were used to produce 19 different solutions in which 
the VLBI antenna positions, source coordinates, clock parameters, and 
environmental parameters were adjusted. For the baselines that were 
determined more than once, a weighted average of the baseline was formed. 

Local survey ties were made by National Geodetic Survey (NGS) at each 
tracking site to tie the VLBI and SLR systems together in a common 
reference system. This survey enables the baselines from one system 
to be calculated with respect to the other system. For these experi- 
ments the survey mark under the laser was the reference point. The 
local survey ties range in length from 16 meters to 8000 meters with 
conservative survey errors being no more than a few parts in 10“6. 

The results of this analysis is summarized in Table 1 where both the 
SLR and VLBI baseline lengths and their differences (SLR minus VLBI) 
are presented. A comparison of the 22 baselines shows a mean diffe- 
rence, of 1.0 ± 1.1 cm with a scatter about zero of 5.2 cm. No 
apparent systematic scale difference betweeen the networks is evident. 

Figure 1 is a map of the baselines and indicates their differences, SLR 
minus VLBI, in centimeters. 


FUTURE EMPHASIS 

Measurements between SLR and VLBI at additional sites in the continental 
U.S. will be made and compared. Of particular interest will be the 
future measurements made for the Monument Peak to Quincy, CA line. SLR 
has over a decade of history measuring across the San Andreas Fault, and 
the results show good agreement with existing tectonic plate theory 
(Minster-Jordan). VLBI has not measured this line frequently enough to 
draw similar conclusions at this time. Measurements for longer base- 
lines, connecting North America and Europe, have been taken and will 
be compared. 
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SLR MINUS VLBI IN CENTIMETERS 



Figure 1. SLR and VLBI Intercomparison Lines 
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GEODYNAMICS OF CRUSTAL DEFORMATION AND SEISMOTECTONIC 
BLOCK MOVEMENTS IN CENTRAL EUROPE 


Han-Shou Liu 


OBJECTIVE 

Geological observations have revealed the style of neotectonic near- 
surface stresses and deformations in central Europe. Seismic 
activity, focal depths and fault-plane solutions of earthquakes 
indicate kinematic reactions within the crust. The objective of 
this study is to develop a crustal deformation which may account 
for the Rhine graben systems and the associated seismotectonic 
block movements in Europe. 

BACKGROUND 

The knowledge of seismic activity and tectonic structure in central 
Europe obtained by the European geoscientists is profound. It 
provides seismic, geological and geophysical facts to test whether 
the solutions of the focal mechanism for earthquakes, the distribu- 
tion of orogenic activity, recent crustal block movements and 
upper mantle structure in this region can be integrated within a 
single Space Geodynamics Program. 

The program for satellite measurement of the Earth's gravity offers 
a variety of models to interpret the interior of the Earth. Among 
the various interpretations, the subcrustal stress field serves as 
an intraplate parameter which corresponds to the lateral variations 
of gravity anomalies. In this study, a computer-aided tomography 
to gravity anomalies is used in determining the crustal stresses in 
central Europe. This study involves a search for tomographical 
interpretations of gravity data with respect to seismic stresses. 

RECENT ACCOMPLISHMENTS 

The subcrustal stresses can be inferred from satellite gravity 
data. Stresses at 6,000 grid points in Europe have been obtained by 
applying tomography to gravity field. Kinematics and dynamics 
have been integrated to show that the measured regional stresses 
in central Europe are derivable from the convection-generated trac- 
tion on the boundary of the elastic spherical shell of the crust 
as inferred from satellite-derived gravity data. The stress pat- 
tern seems to suggest a seismogenic cause for the 1983 Liege earth- 
quake which caused serious economic consequences in Central Europe. 



SIGNIFICANCE 


The stresses as shown in Figure 1 are comprised of several consistent 
stress patterns including the remarkable NW-SE compression under 
central Europe. By- analyzing the sateTTite-deri ved stress fields 
in conjunction with tectonic and seismological observations in 
central Europe, we have obtained the following significant results: 

1. Present-day stress field: 

The present-day stress field in central Europe can be obained 
from focal mechanism of earthquakes, orientation of geological 
joints and in situ stress measurements. The results of these 
seismic and geological observations are summarized in Figure 2. 
Figure 2 shows that the European crust is under NW-SE compres- 
sion which is in full agreement with the satellite-derived 
stress field. 

2. Orientations of the principal compressional and tensional 
stresses and the maximum shear stress in the crust of central 
Europe as inferred from satellite gravity data are shown in 
Figure 3 . These results have provided a dynamical basis for 
the deviatoric stress tensors and displacements in the epi cen- 
tral areas of the historical earthquakes in central Europe. 

3. According to the principles of crystal deformation, the stress 

components in the crustal block of central Europe can be modeled 
as shown in Figure 4. On the basis of this stress field- 

controlled crystal model of crustal deformation, the differing 
seismotectonic regimes of the Upper Rhine graben, the Lower 
Rhine graben and the Belgian and Swabian zones are easy to 
understand. Tensional crustal deformations are occurring 

under the influence of crustal stresses in a distinct manner 
only in the Lower Rhine graben, whereas in the other seismic 
zones, horizontal shear movements are more typical. Figure 4 
shows that seismotectonic block movements in central Europe are 
clearly controlled by the stresses which are inferred from 
satellite-derived gravity data. 

4. This study concludes that it is the stress field as inferred 
from satellite gravity data which determines such phenomena as 
rock stresses in mines, geological joints, earthquake mechanisms 
and seismotectonic block movements in central Europe. The 
crustal stress field in central Europe computed on the basis of 
satellite data gravity resembles closely the experimental 
pattern of single-crystal deformation. It has been shown that 
much of the intraplate deformation in central Europe can be 
explained by the dynamical processes that govern single crystals 
which are a billion times smaller in size than seismotectonic 
bl ocks. 
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FUTURE EMPHASIS 


European seimsologists are largely tutored in certain branches of 
mathematics, and quite prepared for this unorthodox approach to 
seismic activity in Europe. When the European geoscientists could 
not understand the mystery of the Liege earthquake (November 8, 
1983), this approach could probably be used as a last resort. It 
is emphasized that future analysis of stress concentration as 
inferred from satellite gravity data should provide a seismogenic 
model for the strain release pattern in central Europe. 
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Figure 1. Subcrustal Stresses Under Europe as Inferred from Satellite Gravity Data. 
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Figure 2. Present-day Stresses in Centra! Europe Deduced from Seismic and Geological Phenomena. 
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Figure 3. Orientations of the Principal Stresses in the Crust of Central Europe as Inferred from 

Satellite Gravity Data. 
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Figure 4. Seismotectomc Schemeof the Rhine Graben System with the Calculated Crustal Stress 
Components and Possible Horizontal Crustal Block Movements. 
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SENSITIVITY OF SLR BASELINES TO ERRORS 
IN EARTH ORIENTATION 

David E. Smith 
Demos C. Christodoulidis 


OBJECTIVE 

To. assess the sensitivity of inter-station distances derived from 
Satellite Laser Ranging (SLR) to errors in earth orientation. 


BACKGROUND 

The limiting factor in dynamical satellite geodesy for estimating 
station coordinates is clearly our degree of understanding of the 
evolution of the orbit and the precise positioning of the earth 
with respect to the orbit in time. In practice, the laser ranging 
data acquired on Lageos is of sufficient abundance and quality to 
permit the satisfactory simultaneous solution of polar motion, 
earth rotation and station parameters in a least squares adjustment 
procedure (to within an undeterminable absolute longitude 
constant) because each has a unique signal within the laser range 
data. While this procedure has essentially eliminated errors in 
earth orientation as a source of concern, an error propagation has 
been undertaken to fully assess the quality of the results should 
some externally-derived values of polar motion and UT1 parameters 
of unknown error be adopted as a reference within our multi-year 
solutions. 


RECENT ACCOMPLISHMENTS 

An analysis experiment has been performed within our multi-year 
(from the launch of Lageos in 1.97 6 to the end of 1982) solution which 
imposed a known polar motion error on all of the orbital arcs used 
over this interval. All the stations tracking Lageos were permitted 
to adjust and compare to another set of stations derived from a 
set of polar motion data which was assumed to be perfect. The polar 
motion error which was imposed on the recovery contained both long 
and short period errors with a high frequency noi'se-like contribu- 
tion also considered. These errors had an RMS systematic offset 
from the "true" polar motion of about 15 cm with a noise level of 
the order of 10 cm. The errors in UT1 were quite similiar. Since 
this solution was made using the actual contributions from each 
participating station, the experiment included stations which have 
been "fixed" over a good deal of the lifetime of Lageos— our so- 
called "base" stations— with other sites having more limited data 
or transitory occupancies. Therefore, the effect of the averaging 
of these errors over the tracking periods of individual sites was 
also assessed. 
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SIGNIFICANCE 


The SLR technology is capable of relating sites even though they may 
have different tracking intervals with respect to one another. 

These site positions are tied together primarily through the rigor 
of the tracking network comprising the "base" -stations but they also 
require proper earth orientation. When imposing a fixed polar 
motion and earth rotation rate error within this system, the rigor 
(geometric strength) of the network is tested against the misalign- 
ment of the network with respect to the orbital plane used as <a 
reference frame. 

Figure 1 shows the errors in baselines that occurred as a result of 
our perturbing both the pole and universal time from their true 
value by about 15 cm, as described earlier. In general', the num- 
bers range from zero to a maximum of 5 cm and show that for the 
lines connecting our "base stations" the earth orientation errors 
are very small, 1 cm or less. Figure .2 shows .the standard devia- 
tions of the baselines. Generally, these standard deviations are 
comparable to the changes in baseline length shown in Figure 1. 
The larger errors in Figure 1 and the larger s.d.s in Figure 2 
generally correspond to the stations that had a minimum data set. 
For example, station 7896 (Pasadena) was only occupied for effec- 
tively about 6 weeks in late 1980 and during this time the global 
network was weaker than normal. This combination of a relatively 
short stay on site, coupled with a weak supporting network provides 
noticeably weaker baselines that are much more sensitive to error 
in the orbit and models. 

We believe that the following simple relationship between systematic 
errors in earth orientation and the effect of these errors on 
baseline accuracies (o ) summarizes our results: 

Baseline change (error) ~a /-20 x earth orientation error. 

This expression implies that a 10 cm error in earth orientation 
(pole or U.T.) corrupts a strongly determined baseline (a = 2 
cm) by about 1 cm; and weakly determined baseline (a = 5 cm) 
by about 2 to 3 cm. 

The conclusion that we draw from this brief analysis is that base- 
lines between stations that are supported by a global - network of 
tracking stations are only marginally affected by errors in "earth 
orientation. It appears that the global network of stations re- 
tains its integrity even in the presence of systematic changes in 
the coordinate frame, and minimizes the effect of these coordinate 
frame changes on the relative locations of the stations, i.e., 
baselines are not significantly affected. 
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Figure 1. Baseline Error Sensitivity to Earth Orientation. 



1-30 


STN HT 7777777777777777777777777777 
0000000001111111112888889999 
5668899990000111 121333990224 
1232601262359024500345267193 

7051 . 0323322432422323322443342322 QUINCY, CA 

7062 . 302-42 1 2432423223322343442222 OTAY, CA 

7063 . 220321 2322422222222332331221 GREENBELT, MD 

7082 . 3430323433433333333443442432 BEAR LAKE, UT 

7086 . 3223012332422222222343432222 FORTDAVIS,TX 

7090 . 21 12101431311111111222331112 YARRAGADEE, AUS 

7091 . 2223210322422222222333331231 WESTFORD, MA 

7092 . 4434343053433333344443442324 KWAJELEIN 

7096 . 33 23332503423333333332442323 AM. SAMOA 

7102 . 2223212330422222222343332231 GREENBELT, MD 

7103 . 4444434444044444434454444443 GREENBELT, MD 

7105 . 2223212322402222222232331221 GREENBELT, MD 

7109 . 2323212332420222222343341221 QUINCY, CA 

7110 . 3223212332422022222342331221 MONUMENT PEAK, CA 

7112 . 2223212332422202222343331221 PLATTEVI LLE, CO 

7114 , 3323212332422220222343341221 OWENS VALLEY, CA 

7115 . 33232 12332422222022342331221 GOLDSTONE, CA 

7120 , 2223212432322222202332331221 MT. HALEAKALA, HI 

7210 . 22232 12432422222220332331222 MT. HALEAKALA, HI 

7833 . 4334323433423333333055442342 KOOTWIJK, NETH 

7834 , 443 4 423434534444433505542442 WETTZEL, W. GERM. 

7835 . 3323323323423233222550442232 GRASSE, FRANCE 

7892 . 3434433443433333333454043433 VERNAL, UT 

7896 . 4434333443434334333444403333 PASADENA, CA 

7907 . 221221122241 111 11 11222330121 AREQUIPA, PERU 

7921 . 3224212332422222222342431021 MT. HOPKINS, AZ 

7929 . 22232 13223422222222443332202 NATAL, BRAZIL 

7943 . 2212221431311111112222331120 ORRQRAL, AUS. 


Figure 2. Estimate Baseline Uncertainty 
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CHAPTER 2 

GLOBAL EARTH DYNAMICS 


OVERVIEW 


Studies of global earth dynamics are devoted to improving the 
understanding of the origin of motions of the whole earth and 
determining the structural and mechanical properties of the earth's 
interior. The techniques employed in such studies are as varied as 
the studies themselves. They range from determination of changes 
in satellite orbits, to spectrum analyses of polar motion data, to 
studies of the f-lexural strength of the lithosphere. 

1) The first report focuses on the analysis of the earth's rota- 
tion and compares the observations with the atmospheric and 
seismic excitation functions. 

2) The second report continues the study of polar motion. It 
focuses on recent studies which will contribute towards the 
understanding of the excitation mechanism of the Chandler 
wobble. 

3) The third report discusses the use of SEASAT altimeter data to 
obtain information about the flexural strength of the lithos- 
phere. 

4) The fourth report discusses modeling the perturbative force 
due to radiation pressure which acts on the orbit of Lageos. 

5) The fifth report continues the assessment of the perturbative 
force due to radiation pressure which acts on the orbit of 
Lageos. In this latter report however, the diminuation of 
radiation pressure occurring when the Lageos satellite moves 
into the moon's shadow is assessed. 

Contributors to this chapter are B.F. Chao, R. Gross, D.C. McAdoo*, 
and D.P. Rubincam. 


* D.C. McAdoo now is with National Oceanic and Atmospheric 
Administration; National Geodetic Survey. 
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ANALYSIS OF THE EARTH'S VARIABLE ROTATION 


B. Fong Chao 


OBJECTIVE 

To analyze the Earth's rotation, and to compare the observations 
with the atmospheric and seismic excitation functions. 

BACKGROUND 

The variation in the Earth's rotation can be separated into the 
length-of-day (LOD) variation (the axial, or z-component) and the 
polar motion (the x-y component) (Munk & MacDonald, 1960). Possible 
driving mechanisms include: atmospheric/oceanic circulations, 

seismic activities, solar-lunar tides, mantle convection, core- 
mantle coupling, and solar activities (see e.g., Lambeck, 1980).' 
The major problem of concern here is the identification of these 
dynamical processes as primary driving mechanisms for the various 
features in the variation of the Earth's rotation.- 

Besides the solar-lunar tidal terms, the LOD variation can be 
conveniently divided into frequency bands: (1) seasonal and short- 

period (time scale shorter than 1 year, including the dominant 
annual and semi-annual components), (2) interannual (1 to 10 
years), and (3) decade and secular (longer than 10 years). The 
driving mechanism for (1) has recently been confirmed to be the 
changes in the atmospheric angular momentum by Rosen & Sal stein 
(1983), while (2) and (3) have remained largely unstudied. 

The polar motion is composed of the annual wobble, the 14-month 
Chandler wobble and a secular drift. It has been known (Wahr, 
1983) that the energy associated with the atmospheric mass movement 
is sufficient to excite the observed annual wobble, although a time- 
domain correlation study is left to be desired. As for the Chandler 
wobble, despite a great deal of effort by various investigators 
over several decades, the primary driving mechanism is still a 
mystery. The latter is also true for the secular drift. 

The data analysis procedures for the LOD and the polar motion are 
different. This is simply because the Earth's responses for the 
LOD and the polar motion differ in a fundamental way. For LOD, the 
solid Earth acts as an all -pass filter in the sense that whatever 
befalls it, it responds accordingly. This allows us to do a 
direct comparison of the LOD with excitation functions. For the 
polar motion, on the other hand, the Earth is a free oscillating 
system with a natural period of 14 months (the Chandler wobble). 
The transfer function of such a system has a peak at the natural 
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period and a phase delay of 90°. Therefore, prior to any comparison 
of polar motion observations with excitation functions, it is neces- 
sary to deconvolve the observed polar motion— a procedure often 
beset by noise. 

RECENT ACCOMPLISHMENTS AND SIGNIFICANCE 

For LOD, an empirial correlation study was conducted (Chao, 1984a) 
between the LOD variation and the southern oscillation (S0)/E1 Nino, 
the most prominent signal in interannual, global-scale atmospheric 
fluctuation characterized by a seesawing of air mass between eastern 
and western hemisphere in the tropical Pacific-Indian Ocean region. 
Indeed, any atmospheric mass movement, such as SO, that exerts an 
external torque on the solid-Earth will certainly induce LOD varia- 
tions. Thus, two relevant time series for the period 1957-1983 
were compared: the SO Index and the interannual LOD variation. The 
latter is obtained by removing the least-squares estimates for 

the long-period (secular and decade) and- the seasonal variations 
in the BIH (Bureau Internationl de 1'Heure) LOD series. The two 
series have a very encouraging qualitative correlation, in partic- 
ular with respect to El Nino events (especially for the strong 

1982-83 episode); and the linear correlation coefficient is found to 

be 0.55. It is believed that much, if not most, interannual LOD is 
caused by the SO, and the true correlation is considerably higher 

than its apparent value considering the fact that the SO Index is 
merely an indicator derived from two local atmospheric measurements 
while the LOD is a global quantity. 

A comprehensive experimental study of the predictability of the 
polar motion was completed (Chao, 1984b) using a homogeneous BIH 
data set for the period 1967-1983. The numerical prediction of the 
Earth's polar motion is of both theoretical and practical interests. 
Thus, based on a general knowledge of the physics of the annual and 
the Chandler wobbles, a numerical model for the polar motion was 
built by allowing the wobble periods to vary. Using an optimum base 
length of 6 years for prediction, this "floating-period" model, 
equipped with a nonlinear least-squares estimator, is found to 
yield polar motion predictions accurate to within 0"012 to 0"024 
depending on the prediction length up to one year, corresponding to 
a predictability of 91-83%. This represents a considerable improve- 
ment over the conventional fixed-period predictor, which, by its 
nature, does not respond to variations in the apparent -wobble periods 
(in particular, a dramatic decrease in the periods of both the annual 
and the Chandler wobbles after the year 1980). The superiority of 
the floating-period predictor to other predictors based on critically 
different numerical models is also demonstrated. 
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On the theoretical front, a study was conducted (Chao, 1984c) on 
the maximum entropy /'autoregressive (ME/AR-) modeling of time series, 
a. spectral method widely used in geophysics and, in particular, 
the study of the period and damping of the Chandler wobble. Thus,’ 
the AR model of a random process is interpreted in the light of 
the Prony's relation which relates a complex conjugate pair of 
poles of the AR process on the z-plane (or the z-domain) on the 
one hand, to the complex frequency of one complex harmonic function 
in the time domain on the other. It is then seen that the AR 
model of a time series is one that models the time series as a 
linear combination of complex harmonic functions, which include 
pure sinusoids and real exponentials as special cases. An AR model 
is completely determined by its z-domain pole configuration. The 
ME/AR spectrum, defined on the unit circle of the z-plane (or the 
frequency domain), is nothing but a convenient, but ambiguous visual 
representation. The position and shape of a spectral peak is 
determined by the corresponding complex frequency, and the height 
or area of the spectral peak contains little information about the 
amplitude of the complex harmonic functions. 

FUTURE EMPHASIS 

The Southern 'Oscillation (SO)-LOD correlation suggests that the 
atmospheric/oceanic mass transport associated with SO may also 
excite polar motion. In particular, the El Nino events do have a 
time scale close to the natural period of the Earth's wobble— the 
14-month Chandler period. Unless the SO is essentially all 
equatorial, this means an observable signature in the Chandler 
wobble excitation function, which can be derived from the observed 
polar motion time series. A preliminary study using Lageos (LAser 
GEOdynamics Satellite) derived data has shown an encouraging corre- 
lation between the SO Index and the y-component of the polar motion 
excitation function (Gross and Chao, 1984). 
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On 


THE OBSERVED EXCITATION FUNCTION OF THE CHANDLER WOBBLE 

Richard Gross 


OBJECTIVE 

The objective of this research effort is to gain a better under- 
standing of the mechanism that is keeping the Chandler wobble 
excited. 

BACKGROUND 

Any rigid body that is rotating about some axis which is not its 
principal moment of inertia axis will experience a wobble as it 
rotates. The rotation axis will oscillate about the inertia axis 
(as seen in a frame of reference attached to the rotating body) and 
will try to align itself with the inertia axis since their alignment 
represents a state of minimum energy for the rotating body. For the 
earth this motion is known as the Chandler wobble in honor of its 
discoverer. If the earth remained undisturbed during this motion* 
the rotation pole would trace a spiral about the inertia pole 
with a period of 434 days and an e-folding time for the amplitude 
decay of about 40 years. Eventually the motion would damp out 
entirely and the rotation and inertia axes would be aligned. 
However, the Chandler wobble has been under observation for the 
past 80 years and is seen to actually grow in amplitude during 
part of this time. Therefore, there must be some as yet unknown 
mechanism that is exciting the Chandler wobble. Any mechanism 
that causes mass motions of the earth on time scales much less 
than 434 days can be the source of this excitation. Amongst those 
proposed are earthquakes, weather effects and core-mantle 
interactions. 

RECENT ACCOMPLISHMENTS 

The equation relating the position of the rotation pole m to that 
of (essentially) the inertia pole Y is: 


Y(t) = m(t) + — m(t) 


where m(t) is the complex-valued angular displacement of the rota- 
tion pole from the origin of some coordinate system, ¥(t) is 
the similar displacement for the inertia pole, o 0 is the complex- 
valued ei genf requency of the motion and the dot denotes time dif- 
ferentiation. This equation can be solved for m(t) thereby 
obtaining: 
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m(t) = m Q e 


rt icr (t-r) 

icr Q '{'(tje dx+ n(t) 

o 


where m 0 is the initial position of the rotation pole and n(t) 
represents the noise in the observed' time series m(t). By recogniz- 
ing the above integral as a convolution, the problem of recovering 
¥(t) from m(t) i.s seen to be one of deconvolution. In this 
study Y(t) is-- ‘ recovered by performing this deconvolution in 
the frequency domain modeling the noise' as though it were white. 
The result is shown in figures 1 and 2. 

SIGNIFICANCE 

By correlating the features seen in this observed excitation func- 
tion with known physical events on the earth (such as large earth- 
quakes or strong El Nino/Southern Oscillation events) it should be 
possible to determine the source (or sources) of the excitation of 
the Chandler wobble. 

FUTURE EMPHASIS 

Continued interpretation of the excitation function derived from 
various observed Chandler wobble data sets is planned. If each 
and every wiggle in figures 1 and 2 can be physically explained 
then the question of the source of the excitation of the Chandler 
wobble will be solved. 
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SEASAT OBSERVATIONS OF LITHOSPHERIC FLEXURE 
David C. McAdoo 


OBJECTIVE 

In order to understand something about the dynamics of plate 
tectonics it is necessary to know how the lithosphere itself res- 
ponds to tectonic forces. Specifically, we would like to estimate 
the flexural strength of the lithosphere. Further, we wish to 
describe the apparent dependence of this strength on age of lithos- 
phere and magnitude of load. 


BACKGROUND 

Oceanic lithosphere of various ages is subjected to relatively 
large forces as it approaches subduction zones. Regions seaward of 
subduction zones are, in fact, uniquely suited to testing mechanical 
models of oceanic lithosphere. Outer Rises, represent the deforma- 
tional response of the lithosphere to the forces in these regions. 
A substantial number of studies (see McNutt and Menard, 1982 or 
McAdoo and Martin, 1984, for references) have compared bathymetric 
observations of Outer Rises with models of lithospheric flexure. 
This study, which is now in its second year, has tested models of 
lithospheric flexure on SEASAT altimetric observations of the 
geoid over Outer Rises. These altimeter data have been found to 
provide significant new information about the strength of the 
oceanic lithosphere. Among the significant results derived from 
altimeter data is confirmation of the proposition that the effective 
elastic thickness, Tp, of the lithosphere increases with age in 
approximate accord with the relation T s c * age 1 '^. In fact, a 
study (McAdoo and Martin, 1984) of (Juter Rises in six regions 
concludes that altimetric data provide better agreement with the 
age-thickness relation than do bathymetric data (see Figure 1). 
This thickness-age relation can be shown to agree with thermal 
models of the lithosphere when an experimental ly-deri ved ductile 
flow law is invoked for rock comprising the lower lithosphere 
(Goetze and Evans, 1979). In fact, one can interpret effective 
elastic thickness as simply depth to a particular isotherm 
(~700°C) . 
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RECENT ACCOMPLISHMENTS 


The SEASAT study of Outer Rises in six regions — the Aleutian, 
Mariana, Izu-Bonin, Kuril, and Middle American trenches — has been 
accepted for publication (McAdoo and Martin, 1984). This study has 
now been extended to other regions including the Peru-Chile, Puerto 
Rico and Philippine trenches. To date, no bathymetric profiles 
across Outer Rises in these three regions have been analyzed. 
Therefore our altimetric results, shown in Figures 1 and 3, are 
novel. When our results for all nine regions are examined we 
conclude that the lithosphere does continue to strengthen out to 
150 Myr— the approximate age of lithosphere being subducted at the 
Izu-Bonin and Mariana trenches. This conclusion contrasts strongly 
with that of McNutt and Menard (1982) who suggested that the 
lithosphere is weaker than predicted (Goetze and Evans, 1979) from 
experimental deformation of rocks. McNutt and Menard's work was 
based solely on bathymetric data. We find (McAdoo, et al., 1984) 
that the lithospheric rheology, or strength-depth envelopes, pre- 
dicted by Goetze and Evans (1979) and later employed by Bodine et 
al . (1981) are consistent with the altimeter data. Figure 2 shows 
an idealized yield stress envelope for realistic rheology lithos- 
phere (RRL) of two different ages, 60 Myr and 160 Myr. 

Thermally activated creep in the ductile regime governs the shape 
of the envelope (and effective elastic thickness) at depth. A 
nominal strain rate of 10~ 6 s _ l is assumed. Activation energies for 
ductile flow are from Goetze and Evans. Brittle or cataclastic 
failure controls the envelope shape in the shallow regime. An 
elastic core separates the ductile and brittle zones; its thickness 
is load dependent. Using this realistic rheology as sketched in 
Figure 2, a maximum collapse moment can be calculated as a 
function of thickness or age. From the application of the elastic 
lithospheric model to SEASAT data we calculate a zero-crossing 
moment which closely approximates the actual maximum moment (see 
McNutt and Menard, 1982). Therefore, maximum moments from analysis 
of SEASAT data agree well with those predicted by rock mechanics 
experiments (RRL). In summary both the age vs. thickness relation 
and maximum moment vs. thickness relation predicted by rock mecha- 
nics experiments (RRL) agree with SEASAT altimeter observations. 


SIGNIFICANCE 

SEASAT altimeter data over Outer Rises provide an important con- 
straint on mechanical models of the oceanic lithosphere. These 
data are quite consistent with an experimentally predicted mechani- 
cal model (RRL) of the lithosphere which indicates that this model 
may be useful in other geodynamic investigations. 
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FUTURE EMPHASIS 


The potential utility of SEASAT and other satellite altimeter data 
for studying dynamics of the Ti't'hbs'ph'eTe' and sub-lithosphere is 
only beginning to be exploited. 

Similar studies based on SEASAT data will be attempted in intra- 
plate regions. For example, the lithospheric folding in the central 
Indian Ocean Basin can perhaps be explained using the realistic 
rheological (RRL) model of the lithosphere. 


REFERENCES 
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FIGURE CAPTIONS 


Figure 1. Effective elastic thickness, T e , derived from SEASAT 
data (crosses) over nine trenches. Results from Middle America 
(MD), Aleutian (A), Philippine (P), Kuril (K), Mariana (M), and 
Bonin (B), indicated by solid crosses. Results from Puerto Rico 
(PR), Peru (PC) and South Sandwich (S), indicated by dashed crosses, 
T e 's are plotted versus age. Bodine et al . (1981) age-thickness 
relation (9) is shown as curves; see text. Dots are estimates of 
T e derived from ship bathymetric profiles by Caldwell (1978) for 
trenches M, MD, A and P, by Caldwell et al. (1976) for B and by 
McAdoo et al . (1978) for K. 


Figure 2. Yield stress envelope (linearized) for lithosphere (RRL) 
of two different aqes. Adapted from Goetze and Evans (1979). 
Strain rate e = 10“*6 S -1 is assumed (Parsons). 


Figure 3. Maximum moment, M 0 , (actually zero-crossing moment) vs. 
lithospheric thickness as inferred from SEASAT data over 9 Outer 
Rises (see Figure 1). Curve represents prediction from lithospheric 
model (RRL). 
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Figure 1. Effective Elastic Thickness, T e , Derived from SEASAT Data (Crosses) Over Nine Trenches. Results from Middle 
America (MD), Aleutian (A), Philippine (P), Kuril (K), Mariana (M), and Bonin (B), Indicated by Solid Crosses. Results 
from Puerto Rico (PR), Peru (PC) and South Sandwich (S), Incidated by Dashed Crosses, T e 's are Plotted Versus 
Age. Bodine et al. (1981) Age-Thickness Relation (9) is Shown as Curves; See Text. Dots are Estimates of T 
Derived from Ship Bathymetric Profiles by Caldwell (1978) for Trenches M, MD, A and P, by Caldwell et al. (1976) 
for B and by McAdoo et al. (1978) for K. 
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Figure 2. Yield Stress Envelope (Linearized) for Lithosphere (RRL) of Two Different Ages. Adapted from 
Goetze and Evans (1979). Strain Rate e = 10 -16 s -1 is Assumed (Parsons). 
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Figure 3. Maximum Moment, M 0 , (Actually Zero-Crossing Moment) vs. Lithospheric Thickness as Inferred from SEASAT 
Data Over 9 Outer Rises (See Figure 1). Curve Represents Prediction from Lithospheric Mode! (RRL). 
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LAGEOS ORBIT AND THE ALBEDO PROBLEM 
David P. Rubincam 


OBJECTIVE 

To obtain an analytic expression for the radiation pressure force 
on a satellite due to sunlight reflected from the earth. (This is 
commonly called the albedo problem). 

BACKGROUND 

The Lageos satellite undergoes unexplained along-track accelerations 
on the order of 3 x lO" 12 ms** 2 , once the average charged particle 
drag effect is removed. These accelerations are believed to be due 
mainly to terrestrial radiation pressure (Smith, 1983; Anselmo, et 
al . , 1983). The effect of sunlight, reflected off the surface of 
the earth must thus be modeled to insure an accurate orbit for 
Lageos. An accurate orbit is necessary for carrying out Lageos' 
mission of measuring tectonic plate motion, polar motion, and earth 
rotation. 

An analytic expression for the force due to reflected sunlight is 
desirable to avoid a time-consuming numerical integration when 
computing the orbit. Obtaining such an analytic expression is a 
venerable problem in celestial mechanics (e.g., Lockry, 1966; Smith, 
1970, and Lautman, 1977a, 1977b). Difficulties include integrating 
over only the sunlit portion of the earth seen by the satellite, 
assuming that the surface follows Lambert's law, and a geographi- 
cally and temporally varying albedo. 

RECENT ACCOMPLISHMENTS 

The present investigation focuses on a spherical harmonic approach 
to the problem. An equation for the force has been obtained by 
assuming the earth's surface reflects sunlight according to 
Lambert's law. The equation is an integral over the whole earth's 
surface. Expressions occurring inside the integral are expressed 
in terms of spherical harmonics. The problem is thus reduced to 
integrating products of spherical harmonics. 

FUTURE EMPHASIS 

Future work will concentrate on working out the integration of the 
product of the spherical harmonics. The integration will give a sum 
of a large number of terms. The dominant terms for Lageos will be 
found and their effect on the orbit will be assessed. 
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LAGEOS ORBIT AND SOLAR ECLIPSES 


David P. Rubincam 


OBJECTIVE 

To assess the importance of solar eclipses on Lageos ' orbit. 
BACKGROUND 

Solar radiation pressure perturbs the orbit of the Lageos satellite. 
The GEODYN orbit determination computer program includes solar 
radiation pressure as one of the forces operating on the satellite 
as it integrates the orbit. GEODYN also takes into account the 
extinction of sunlight when Lageos moves into the earth's shadow. 

GEODYN does not at present take into account the diminution of 
radiation pressure when Lageos moves into the moon's shadow, i.e., 
suffers a solar eclipse by entering the moon's umbra or penumbra. 
This diminution will affect Lageos' orbit by weakening for a time 
the radiation presure acting on the satellite. The importance of 
this effect must be assessed in order to assure an accurate integra- 
tion of Lageos' orbit. An accurate orbit is necessary for Lageos 
to accomplish its mission of monitoring tectonic plate motion, 
polar motion, and earth rotation. 

RECENT ACCOMPLISHMENTS AND SIGNIFICANCE 

The effect of solar eclipses on the semimajor axis of Lageos' orbit 
has been computed analytically by assuming Lageos to be in a circular 
orbit, the sun and the moon to be in the plane of the orbit, and 
the moon to be stationary in the sky in front of the sun. Also, 
the magnitude of the radiation pressure is assumed to be linearly 
related to the angular separation of the sun and moon, and that 
Lageos is a perfect absorber of radiation. 

The computation indicates that an eclipse of the sun by the moon as 
seen by Lageos can affect the semimajor axis at the 1 centimeter (1 
cm) level. {On the time-scale of a day this will appear as an 
instantaneous change in the semimajor axis — an eclipse will last 
about an hour at most). Such a change is significant enough to 
include in GEODYN, in order to get an accurate orbit for Lageos. 

FUTURE EMPHASIS 

The task ahead is to determine just how often Lageos has been 
eclipsed in the past and how each eclipse affected the orbit. To 
do this, a computer program will be developed to look at the solar. 
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lunar, and Lageos positions around the time of each new moon since 
Lageos' launch. It will see if an eclipse took place, and if s.o.,. 
the perturbation of the orbit. About 20. penumbral and no umbral 
eclipses are expected. 
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CHAPTER 3 

GRAVITY FIELD MODEL DEVELOPMENT 


OVERVIEW 


Knowledge of the earth's gravity field is fundamental to under- 
standing the dynamics of the earth. For solid earth physics, 
knowledge of the variations in-the gravity field provides informa- 
tion on the earth's physical properties and geodynamic processes 
and places constraints on the internal structure of the earth. In 
oceanography, knowledge of departure of the actual sea surface 
from a unique equipotential surface of the earth's gravity field 
(the geoid) can reveal information on oceanic circulation. In 
addition, other areas which benefit from knowledge of the earth's 
gravity field are satellite orbit determination and classical 
geodesy. This chapter describes: 

1) The development and evaluation of a specialized gravity field 
model for Lageos which is used for the analysis of satellite 
laser ranging measurements and the accurate relative changes 
in these locations are used to study determination of tracking 
station locations of tectonic plate motion. 

2) The developmental status of the Geodynamics Computational 
System (GEODYN) 

and 

3) The development of ancillary computer software supporting 
geodynamic research in terms of error analysis and eigen-vector 
analysis. 

Contributors to this chapter are F.J. Lerch, B.H. Putney, and 
T.L.- Felsentreger. 
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VERIFICATION OF THE ACCURACY OF GEM-L2 
IN RESPONSE TO CRITICISM BY LAMBECK AND COLEMAN 

Francis J. Lerch 


OBJECTIVE 

The objective is to evaluate the accuracy of GEM-L2 and its 
improvement in estimating Lageos orbits in order to provide better 
baselines for plate tectonics, improved polar motion and earth 
rotation. Analyses and comparisons with other models are made to 
verify our accuracy in contrast to the Lambeck and Coleman result 
which denied our accuracy, estimates. Also the purpose of GEM-L2 
is to improve the low degree and order geopotential which will 
provide for refinement of our knowledge of the broad features of 
ocean circulation. 


BACKGROUND 

The development of the gravity field for Lageos was a principal 
investigation in the Crustal Dynamics Project. Based upon the 
previous Goddard Earth Model, GEM 9, a global set of baselines 
(Table 1) were improved from 7 cm with this model to 1.8 cm with 
the. new GEM-L2 model. However, the accuracy (Table 2) of GEM-L2 
was claimed in a recent publication by Lambeck and Coleman (1983) 
to be no better for the low degree terms than the differences with 
other existing models (GRIM 3, SAO 77). These differences, particu- 
larly degree 3 are an order of magnitude larger (Table 3) than the 
accuracy of GEM-L2. Evaluation and comparisons are made to dispute 
the Lambeck /Coleman claim and to show that their "error spectrum" 
estimate for models is. highly unreliable. 


RECENT ACCOMPLISHMENTS 

A refined gravity field model, Goddard Earth Model GEM-L2 (Lerch, 
et al. ), has been derived using the Lageos orbital data yielding 
better baseline measurements (Table 1) for the analysis of tectonic 
plate motion. This field also contributes to a better under- 
standing of the broad features of ocean circulation through its 
improvement of the long wavelength geoid (terms through -degree and 
order 4) with an accuracy estimated at ± 8 cm (Table 2). In GEM-L2 
two and half years of Lageos laser data acquired from over 20 
well -distributed stations were combined with the existing data 
from the best satellite-derived model, GEM-9 (Lerch, et al.). The 
accuracy estimates of GEM-L2 are verified below. 
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First the general accuracy of the model is reviewed from an evalua- 
tion based upon new gravity anomaly data (Table 4). Letting the 
symbol <x> denote an average value of x over the earth Kaula's 
formula (1966), namely 

EST cr 2 D/\y/\(Ag s ) = <Ag s > - <Ag s Ag> 

was used as a test criteria to estimate the accuracy (commission 
error) of the gravity anomaly (a s ) of a satellite model from a 
global distribution of gravity anomaly data (a ). The gravity 
anomaly error in GEM-L2 may be computed from the coefficient errors 
as follows: 

02GEM-L2 U9S) ' T * 2 I (n-DH^tV + oMSJ) 

n=2 fifO m 


This error estimate is then compared in Table 4 with Kaula's test 
criteria above to verify the GEM-L2 error estimates using a number 
of gravity anomaly data sets. A calibration factor k is defined 
in the table and was employed to compare the agreement of the 
error estimates, namely 

EST „D AT A < A 9 s > “ ko GEM-L2 (4g s> 

The satellite a has already been calibrated upward by a factor 
of /10 as previously done since the GEM-5 model (Lerch, et al., 
1974) using the same method with surface gravity data only. The 
new gravity anomalies, using SEASAT over ocean areas and gravi- 
metry over land areas, show that the errors in the satellite model 
are verified and in fact should be reduced by about one-third (k - 
.67) from the previous estimate which was determined from the 
worldwide use of gravimetry only (k = 1.0). The improved set of 
gravity anomalies is responsible for the better calibration factor. 

Since the gravity anomaly data are more sensitive to the higher 
degree terms, a test was made to estimate the calibration factor 
for the errors in the low degree terms to which Lambeck and Coleman 
primarily objected. In Table 5 independent solutions from 
optical data only and laser data only were used to calibrate the 
low degree error estimates for terms out through degree and order 
10. As described in the table the coefficient differences are 
compared with satellite cr's for each degree H. The result shows 
the same agreement from the calibration factor (k = .67) as obtained 
from the new gravity anomaly data above in Table 4. 

In Table 3 a comparison of the 3rd degree coefficient differences 
between GEM-L2 and the fields of SAO 77 and GRIM 3 verify the 
Lambeck/Coleman statement that these are an order of magnitude 
larger than the accuracy estiamtes of GEM-L2 in Table 2. However, 
a more recent European companion field to GRIM-3, namely GRIM-3B, 
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has an rms difference of 10 cm with GEM-L2 which is similar to GEM 
9 as in Table 2. This result with GRIM-3B nullifies the argument 
that recent models from different institutions do not agree closely 
with one another. 

A comparison f.or low- degree -terms between data on 24-hour satellites 
and different gravity models was obtained by Carl Wagner (1983) and 
presented in Table 6. These results support the accuracy of t'he 
GEM-L2 model and show that the large differences with the SAO and 
GRIM-3 models found above correspond to errors in these models 
rather than to GEM-L2. 

It is now shown that the error spectra estimator of Lambeck 
and Coleman for terms of degree £ is highly unreliable. This is 
particularly so for low degree terms which they used to evaluate 
our fields (Tables 1,2, 4, 5, and 6 of their report). The estimator 
is given by (equation 8C of their report) 



where Cn is the field to be tested and C a is a reference field,. and 
where tne vector C £ denotes the spherical harmonic coefficients of 
degree £ and C^ the power in these terms. Using a to denote the 
power of the coefficient errors of degree £ then the expected 
value of e^ 2 (assuming the errors in C& and Tjj, are unbiased and 
uncorrelated) is 

£ 

E K> ■ ‘i ■ 2 r° 2 < c »> + ° 2 

Hence the estimator is unbiased. 

The variance of the estimator represents the spread or deviation 
of the estimator from its expected value, namely 

tf 2 (e£) = E[e 2 - E(e 2 )] 2 

I£ is used to measure the reliability of the estimator ef, i.e., if 
a (e 2 ) » e| the estimator is significantly unreliable and if 

a 2 ( e 2 )< e a it is considered a good estimator. 

£ £ 

A simplified case can be used to test the estimator. Let individual 
errors in the coefficients of degree £ in the test field be normally 
distributed with equal variances (a ? =cx|/2jf+l) and let the power in 
individual coefficients be given by Kaula’s rule (k 2 = 10 -10 /£ 4 )* 
Then the predicted error spread is given by £ 



2£+l 
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OF POOR QUALITY 


_ 2 

where a ^ is the error power spectra for the reference field. 
Denoting the linear spread in the estimator as 

c(e,) = MeJ)] 4 

a reliability factor for the estimator is given as follows 


where r denotes the root square error rati o„ of the variance of 
the estimator (e^ ) to the expected value (a^) of the estimator. 

Two cases are given in Table 7 to show the reliability r of the 
estimator: 


r = — 
o cr 




2£+l- 


for = 0 



f°r = a t 


For purposes of comparison with results from Lambeck and Coleman 
given in Table 6, the error spectra are taken from GEM-L2. 

From inspection of the table this simplified statistical case for 
r=r 2 shows that the estimator significantly over estimates the 
errors in GEM-L2 by a factor of r = 60 for t=2 to r = 9 for 
£=5 and the predicted spread a(fe«) corresponds reasonably with the 
Lambeck /Coleman estimator. Hence the estimator essentially has no 
reliability for estimating the error power spectra of GEM-L2. Also 
note that even if the reference field is perfect (u£ = 0) the 
estimator has no reliability in the error spectra for Jt = 2 to 5 
corresponding to the low degree terms. 

The rms of laser residuals on Lageos arcs are presented in Table 8. 
The residuals of the SAO 77 field are an order of magnitude larger 
than GEM-L2 values. This serves to verify that the differences in 
the coefficients between GEM-L2 and SAO 77, which are an order of 

magnitude larger than the errors in GEM-L2 (see Tables 2 and 3), 

correspond to errors in the SAO field and not to those in GEM-L2 
as claimed by Lambeck and Coleman. It is of further interest to 
note in Table 8 that GEM 9, which is independent of the Lageos 
data, does just as well with the Lageos orbital residuals as 

GRIM-3B and the University of Texas fields which employed this 

data in their solutions. 
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SIGNIFICANCE 


The analysis presented here completely verifies the accuracy of the 
GEM-L2 model and disproves the statistical methods of Lambeck and 
Coleman rejecting the accuracy of GEM-L2 as published in their 
report. The baselines derived from Lageos with GEM-L2 have been 
estimated to have a 2 cm uncertainty due to the errors in the GEM- 
12 gravity model which is supported by the results in Table 1. The 
calibration tests of this report indicate that the errors in GEM-L2 
should be reduced by about 30 percent which is in the opposite 
direction to the result obtained by Lambeck and Coleman. 


FUTURE EMPHASIS 

Future work is continuing to improve the Crustal Dynamcis Gravity 
Model, GEM-L2, with more recent and accurate laser data. More 
extensive analyses verifying the accuracy of the GEM-L2 model will 
be published in a Special Lageos Issue of JGR. 
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Table 1 

Baseline Comparison from Two Independent Lageos Data Sets 
Spanning 1979 Through 1980* 



GEM-9 



GEM-L2 


w/BIH** 

w/LAGEOS 

w/BIH 

w/LAGEOS 

w/LAGEOS 


Polar Motion 

Polar Motion 

Polar Motion 

Polar Motion 

Polar Motion 


and A1-UT1 

and A1-UT1 

and A1-UT1 

& BIH A1-UT1 

and A1-UT1 

Baseline Agreement 




✓ 


for Eight "Base" 
Stations 

{All 28 Baselines) 

7.2 cm 

6.0cm 

» J 

6.9 cm 

4.7 cm 

1.8cm 

"Base" Stations are: 

7063 GSFC, MD 
7086 Ft. Davis, Texas 

7090 Yarragadee, Australia 

7091 Westford, MA 
7115 Goldstone, CA 
7114 Owens Valley, CA 
7907 Arequipa, Peru 
7943 Orroral, Australia 


*Thedata has been divided into alternating 15 day segments and therefore both solutions have similar and uniform data distribution from 
these sites. 

**BIH circular D 90-day smoothed values. 












Table 2 

Estimated Individual Low Degree and Order Coefficient Accuracy (cm) 


ACTUAL COEFFICIENT 



UNCERTAINTY 

UNCERTAINTY 

DIFFERENCE (GEM-9) 

COEFFICIENT 

GEM-9 

GEM-L2 

MINUS (GEM- 

-L2) 

(2, m) 

ffJ £, m 


AC 

AS 

C2,0 

0.4 

0.3 

- .3 


C3.0 

1.0 

0.1 

.5 


C4,0 

0.8 

0.7 

.2 


C,S2,1 

2.3 

1.3 

.5 

.6 

2,2 

3.1 

0.9 

- 2.5 

.9 

3,1 

4.8 

2.3 

- .4 

1 .- 

3,2 

7.2 

3.4 

- 7.7 

- 4.1 

3,3 

8.6 

3.2 

-13.1 

- 2.4 

4,1 

4.5 

3.7 

1.5 

.8 

4,2 

4.4 

3.3 

- 1.3 

1.2 

4,3 

4.1 

2.4 

- 3.0 

- .3 

4,4 

5.6 

2.4 

- 2.9 

- 5.1 

TOTAL: RSS GEOID ACCURACY FOR 




4x4 

15.9 cm 

8.1 cm 

17.7 cm 



GLOBAL 4x4 GEOID 2800 CM 
PERCENT ACCURACY OF GEM-L2 = 0.3% 



Table 3 

Comparison of 3rd Degree Coefficients of Recent Gravity Fields with GEM-L2 



COEFF 

VALUES 

GEM-L2 

x 10 6 

SAO 72 

- COEFFK 
SAO 77 

2IENT DIFFERENCES (CM) - 
(GEM-L2- OTHER) 

GRIM 2 GRIM 3. 

GRIM 3B 

C3,0 

0.958 

- 1.5 

- 1.5 

- 1.9 

- 0.6 

-0.0 

C3,1 

2.029 

20.0 

-12.8 

42.7 

- 5.1 

1.9 

S 

0.250 

17.0 

-17.2 

60.5 

1.3 

1.3 

C3,2 

0.903 

79.5 

- 9.3 

24.8 

-50.6 

-7.8 

S 

-0.616 

88.8 

41.5 

-45.3 

-57.8 

-3.1 

C3;3 

0.722 

148.0 

36.2 

6.4 

13.4 

-5.1 

S 

1.414 

-72.7 

-47.2 

-91.2 

-52.3 

-3.1 

RSS DIFFERENCE (CM 

205.1 

76.2 

128.5 

94.1 

10.6 
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Table 4 

Calibration of Satellite IVlodel Errors 


EST a Data ^g s ) - k ^Model {Ag s) = Formal (Ag s> 

TEST CALIBRATION FACTOR K =yf\Q IN GEM L2 (Also used in GEM 5-9) 

ctqEM l2 ( a 9$) = A -63 m 9als Commission Error 
GRAVITY ANOMALY SATELLITE MODEL ERRORS CALIBRATION 


DATA TEST CASES (Ag*) • 

EST °Data (Ag s) 

FACTOR 

5° Equal Area Surf. Grav. (1976) 

5.1 mgals 

1.1 

5° Equal Area Surf. Grav. (1981) 

4.8 

1.0 

Ag ALTIMETER (OCEAN) 



1° x 1° SEASAT Altim. + S.G. Land (1981) 

3.1 

0.67 

5° x 5° SEASAT Altim. + S.G. Land (1981) 

3.0 

0.65 

5° x 5° SEASAT Altim. Only (1981) 

3.0 

0.65 

5° x 5° GEOS-3 Altim. Only 

3.2 

0.69 


*Data Sets Obtained from Rapp. 
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Table 5 

Calibration of Satellite Model Error Estimates from Independent Solutions of 
Sat Opical Data Only and Laser Data Only 


CALIBRATION (n = degree 2) 

RMS^ * (OPT.) + a 2 (LASER)] 

RMSp = RMS of difference in coefficients of degree n between OPT. 
and laser solutions (n = fi) 

k n = Calibration factor 

n 

a n ~ ^(C^) + (^nm^ /2n+1) 

m=0 

Satellite Coefficient Errors 
cr(C,S) -y/10 CT PORMAL (C,S) 


UNITS OF 10~ 9 


DEG n 

RMS n 

ff n (OPT) 

a n (LASER) 


3 

14 

18 

11 

.66 

4 

10 

10 

13 

.62 

5 

17 

21 

23 

.55 

6 

15 

15 

18 

.66 

7 

24 

25 

24 

.69 

8 

19 

19 

26 ‘ 

.59 

9 

24 

25 

30 

.62 


AVERAGE .63 
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Table 6 

Gravity Model Comparison Using Independent Data from 24 Hour Satellites 


EUROPE 


Weighted J5E1 SE2 SE4 SE6 GRHVI3 GRIM3B GEM3 GEIVI7 GEM9 GEIVUOB GEML2 
RMS 

Residual i [ 

Accelerations 63.4 18.7 51.7 18.0 25.7 4.2 . 7.4 4.1 5.9 4.2 1.0 


where the weighted residual acceleration is (A 0 — A c /cf^) 

Contribution of Geopotential by Degree to Acceleration 


Contribution 
from Full 
Field 


Contribution 
from Estimated 
Error of GEM-L2 


3000x 10* 
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Table 7 

Reliability of Lambeck/Coleman Error Estimate 


RELIABILITY FACTOR 



r o 

r 2 

DEGREE £ 

ii 

o 

a £ = 2(7 £ 

2 

40 

60 

3 

12 

18 

4 

3 

12 

5 

6 

9 

6 

3.3 

5 

10 

1.3 

2 


PREDICTED 

SPREAD 

LAIV1BECK/ 

COLEMAN 

.ERROR 

ESTIMATE 

ERROR 

SPECTRUM 

17 eg = r 2 CT £ 

eg(GEM-L2)* 

agGEM-L2) 

±120 

149 cm 

2 cm 

± 90 

-197 cm 

5 cm 

± 72 

— 91 cm 

6 cm 

±126 

187 cm 

14 cm 

± 80 

16 cm 


±' 70 


35 cm 


*Table 4 results from Lambeck/Coleman report. Reference Field SAO 1980. 
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Table 8 

Gravity Model Comparison 


LAGEOS LASER RANGE DATA 
RMS OF FIT - THREE 7 DAY ARCS 
(CM) 



2/1-9 

2/9-16 

2/22-30 

TOTAL 

GEM-L2 

10 cm 

12 

8 

10 

GRIM3B 

23 

19 

22 

21 

UNIV. OF TEXAS 

25 

22 

20 

22 

GEM 9 * 

21 

19 

20 

20 

SAO 77 * 

113 

105 

128 

114 


* GRAVITY FIELD INDEPENDENT OF LAGEOS DATA. 
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GEODYN SYSTEMS DEVELOPMENT 


Barbara H. Putney 


OBJECTIVES 

The purpose of the Geodyn Orbit Determination and Parameter Estima- 
tion, the SOLVE and ERODYN Programs is to recover geodetic and 
geophysical parameters from satellite and other data in a state- 
of-the-art manner. 

BACKGROUND 

In 1971 the NONAME and GEOSTAR programs were combined to create the 
GEODYN program. The SOLVE program was created at the same time. 
A few years later the ERODYN, error analysis program was written. 
The philosophy of the development of the software system has been 
maintenance of computer efficient, well -structured software, with 
appropriate orbit, earth and numerical models, using precise satel- 
lite measurement modeling and efficient numerical models, and 
performing careful benchmark procedures. This care has paid off in 
the production of several GEM's (Goddard Earth Models), precision 
station locations, improved tidal, GM, polar motion and earth 
rotation values, consistent baselines, and encouraging GRM simula- 
tions. Careful usage, analysis, and modeling using laser, alti- 
meter and other satellite data from LAGEOS, SEASAT, STARLETTE, 
GEOS and BE-C satellites as well as many others has made these 
accomplishments possible. 

RECENT ACCOMPLISHMENTS AND SIGNIFICANCE 

Continued solutions for gravity field, pole positions, earth rota- 
tion, GM, and baselines have been made as part of the Crustal 
Dynamics Project. Some tidal parameters have been recovered as 
well. The eight digit station identification number has been 
incorporated in the software and new techniques for constraining 
monthly station parameters to each other is being developed. This 
is allowing the analysts even more flexibility in the shaping of 
solutions from monthly sets of normal equations and right-hand 
sides. 

The GRM simulations have satisfied the analysts that with the 
simplistic model of the gravity field, GRM accurate gravity models 
can be obtained. It is expected that the Cyber 205 will enable 
finer detailed gravity models to be used. 



Conversion of the Geodyn and Solve programs to the Cyber 205 computer 
has continued this year. The Geodyn II (vector rewrite) program 

had its first successful orbit generator run on Oecember 13, 
1983. Comparison running times are included, in -Figure 1. The' 
initial timings ar.e very encouraging. The next task is to process 
a data arc. The program is designed to run on both the IBM/Amdahl 
or the Cyber 205 computer, sharing a common source file and 

requiring a utility program to create the individual source 
programs. 

This year the Solve program inversion routines have been fully 
vectorized and that portion of the program is running 25 times 
faster than the IBM 3081. It is truly remarkable that the Cyber 
205 can fully invert the 1921 x 1921 matrix in 7 minutes or less. 
The details are in Figure 1. There is an I/O problem with this 

program and currently we are achieving only 33 percent CPU utli- 

zation. Two things can improve this situation, increased real 
memory (currently our Cyber only has one million words of memory) 
and concurrent I/O usage. We are working on the second item* but 
with only 2 channels on this system the best one can hope for is 
50 percent CPU utilization. It is hoped that additional memory 
will be purchased for the machine. Work continues on vectorizing 
the data management areas of the program. 

In addition to the Geodynamics/Crustal Dynamics Programs, the 
Geodyn software is being used by the AMPTE Project. In support of 
this effort this year, the Geodyn program has been sent to the 
Applied Physics Labortory in Laurel, Maryland; DFVLR in Germany; 
JPL in Pasadena, California. In support of Geodynamics/Crustal 
Dynamics research the program has been sent to Australia, The 
Netherlands, Japan, and Ohio State University, 

FUTURE EMPHASIS 

The Topex gravity field will require complete check out of existing 
models and the addition of others, such as sea surface topography. 
This is going to be a large scale data processing effort and the 
software will be made as efficient as possible for this purpose. 

The Geodyn II program needs to be completed and debugged for all 
currently implemented data types. New data types such as GPS should 
be implemented as well. 

The Solve program needs to become a more flexible tool for Crustal 
Dynamics users. The handling of baselines, station solutions, polar 
motion and earth rotation solutions needs to be presented in a better 
way to make it easier for the analyst to easily distinguish 
between obviously poor and more promising solutions. 
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The Geodyn II program will be modified to have the year 2000 as the 
reference coordinate system and to achieve the other MERIT 
standards. 

Every effort will be made, as resources allow, to continue to mold 
the software in a timely manner to support the needs of the scien- 
tific community determining earth models. 
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GEODYN TEST ORBIT GENERATOR RUN 
2 DAY, 60 SEC STEP SIZE 

36 x 36 GRAVITY MODEL, THIRD BODY EFFECTS, 
STANDARD TIDES AND SOLAR RADIATION 


AMDAHL V-7 

GEODYN I 55 CPU SEC 

GEODYN IIS 1 .8 CPU SEC 

GEODYN HE 80 CPU SEC 

GEODYN HE 

(ALL SMALL PAGES) 

, CYBER 7-10 TIMES FASTER THAN V-7 

SOLVE TEST-FULL INVERSION 

1921 PARAMETERS - 1 MATRIX INPUT 


IBM 3081 


SOLVE 8201.1 (9/82) 
SOLVE 8308.0 (9/83) 


116.07 CPU MIN 
31.45 I/O MIN 


CYBER 25+ TIMES FASTER THAN 3081 


Figure 1. Time Comparisons of Vector/Scalar Geodyn and Solve Programs 


CYBER 205 


5.4 CPU SEC 
7.9 STU SEC 


CYBER 205 


122.9 CPU SEC 
356.1 STU SEC 


as of 2/84 


3-18 



V 



30468 



THE ERODYN AND QRPIG COMPUTER PROGRAMS 
Theodore L. Felsentreger 


OBJECTIVES 

The role of the ERODYN computer program is to provide error analyses 
involving orbital, geodetic, and geophysical parameters. It has 
been designed to operate as a companion program to the GEODYN orbit 
determination and parameter estimation program. The j) R_ Partitioned 
Ei genvalue/Eigenvector analysis program (QRPIG) is designed to 
process symmetric matrices with an out-of-core partitioning algo- 
rithm for the eigenvectors and eigenvalues of the matrices. 


BACKGROUND 

The initial development of ERODYN was begun in the early 1970s, a 
few years after the GEODYN and SOLVE programs were created, and'was 
intended to be a companion program to GEODYN. Basically, the 
purpose was to provide information on the effect of (assumed or 
actual) errors in certain parameters (unadjusted) upon the values 
of other parameters (adjusted) to be solved for. To accomplish 
this, ERODYN requires as input the matrix of normal equations ("E"- 
matrix) and, optionally, the matrix of variational parti als (“V" 
matrix) as computed within GEODYN. In order to examine various 
partitionings of the total parameter set into adjusted and unad- 
justed parameters, the initial GEODYN run must have all parameters 
of interest "adjusted"; however, only one estimator iteration is 
required using "nominal" initial conditions to generate the required 
normal matrix and variational partials. These matrices are then 
subdivided within ERODYN into the appropriate adjusted and unad- 
justed partitions as specified by the ERODYN user. 

As GEODYN 1 s capabilities have been expanded and improved, an effort 
has been made to maintain the compatibility of ERODYN with GEODYN. 
In addition, many improvements have been made resulting in much 
more efficient use of computer resources. One of these improvements 
has been the application of a different algorithm for the 
eigenvalue/eigenvector option,, resulting in more speed and 
reliability. A product of this change is a stand-alone eigenvalue/ 
eigenvector program called QRPIG. 
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RECENT ACCOMPLISHMENTS AND SIGNIFICANCE 

The GEODYN/ERODYN error analysis approach is extremely efficient in 
terms of computer resource utilization. A typical analysis problem 
would require multiple runs of an error analysis program using the 
same tracking schedule and problem geometry, but with different 

partitionings of the parameter set into adjusted and unadjusted 
parameters. However, with GEODYN/ERODYN the operations which are 
common to these multiple runs are performed only once by GEODYN and 
stored. ERODYN then does the matrix partitioning and linear algebra 
required for each individual error analysis run. In some instances, 
the "RSS iteration" feature of ERODYN can be used to achieve in a 
single run results which would usually require several. 

The GEODYN/ERODYN error analysis approach also minimizes software 
design requirements. Any new measurement type or force model 
modification added to GEODYN is automatically included in the 
error analysis. This eliminates the requirement for separate 
coding, checkout and validation of both the error analysis and 
orbit determination programs. 

The previous method for computing eigenvalues and eigenvectors 
within ERODYN employed the Jacobi algorithm and required that 

sufficient core storage be available to hold the entire matrix. 

The newly-implemented method uses an out-of-core partitioning 
algorithm, first reducing the (real, symmetric) matrix to a 
tri diagonal symmetric matrix using Householder transformations 
(Martin, et al , 1968), and then employing an iterative QR method to 
obtain the eigenvalues and eigenvectors (Bowdler, et al , 1968). 

This whole routine has been packaged as a stand-alone program 
(QRPIG). 


FUTURE EMPHASIS 

An effort is now underway to convert the ERODYN program to the Cyber 
205 machine. The next step will be to optimize the program by 
making efficient use of the vector capabilities of the machine. It 
is anticipated that QRPIG will undergo the same conversion. 

ERODYN can be an important tool for such missions as TOPEX and GRM, 
which will be quite sensitive to the gravity field of the Earth. 
Sensitivity analyses succh as those performed by ERODYN are extre- 
mely valuable in any research involving determination of Earth 
model parameters. 
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